## MJ
## 2017

## Start Matter ------------------------------------------------------

## Load packages.
library(sandwich)
library(lmtest)
library(broom)
library(ggplot2)
library(gridExtra)

## Load functions.
source("polcom-conjoint-utils.R")

## Example 1: Trust conjoint -----------------------------------------

## Get data.
trust <- read.csv("data-trust.csv", stringsAsFactors = FALSE)

## Fix spelling error.
trust$source_adds <- gsub("Have", "Has", trust$source_adds)
trust$source_comments <- gsub("Have", "Has", trust$source_comments)
trust$source_status <- gsub("newpaper", "newspaper", trust$source_status)

## Labels to use for treatment headers in the figure.
trust_labels <- data.frame(
  treatment = c(
    "source_adds",
    "source_party",
    "source_ethics",
    "source_comments",
    "source_readership",
    "source_status",
    "source_soft",
    "source_established"
  ),
  treatment_label = c(
    "Advertisement",
    "Editor/owner openly supports party",
    "Violations of press ethics last year",
    "Comments field",
    "Readership",
    "Current publication form",
    "Entertainment news",
    "When founded"
  ),
  stringsAsFactors = FALSE
)

## Make the vars factors.
for (var in trust_labels$treatment)
  trust[[var]] <- factor(trust[[var]])

## Fix treatment levels and order as per wanted in the figure.
trust$source_adds <- factor(trust$source_adds, levels = levels(trust$source_adds)[c(3, 2, 1)])
trust$source_comments <- factor(trust$source_comments, levels = levels(trust$source_comments)[c(3, 2, 1)])
trust$source_established <- factor(trust$source_established, levels = levels(trust$source_established)[c(1, 3, 2)])
trust$source_ethics <- factor(trust$source_ethics, levels = levels(trust$source_ethics)[c(3, 1, 2)])
trust$source_party <- factor(trust$source_party, levels = levels(trust$source_party)[c(6, 8, 10, 9, 1, 4, 5, 2, 3, 7)])
trust$source_readership <- factor(trust$source_readership, levels = levels(trust$source_readership))
trust$source_soft <- factor(trust$source_soft, levels = levels(trust$source_soft)[c(1, 3, 2)])
trust$source_status <- factor(trust$source_status, levels = levels(trust$source_status)[c(2, 3, 1)])

## Estimate the AMCEs.
trust_amce <- amce(
  data = trust,
  post = "source_choice",
  treatments = c(
    "source_adds",
    "source_comments",
    "source_established",
    "source_ethics",
    "source_party",
    "source_readership",
    "source_soft",
    "source_status"
  ),
  cluster = "resp_id"
)

## Adding a row with treatment labels so the plot have a header on the y axis.
trust_amce <- trust_amce[rev(order(trust_amce$value)), ]
trust_amce$value <- paste0("plain('", as.character(trust_amce$value), "')")
trust_amce$value <- factor(trust_amce$value, levels = unique(trust_amce$value))
add_data <- trust_labels
names(add_data)[2] <- "value"
add_data$value <- gsub("\n", " ", as.character(add_data$value))
add_data$value <- paste0("bold('", add_data$value, "')")
for (var in c("value_order", "estimate", "std_error", "statistic", "p_value")) {
  add_data[[var]] <- NA
}
add_data$value_order <- -Inf
trust_amce <- rbind(trust_amce, add_data)
trust_amce <- trust_amce[order(trust_amce$treatment, trust_amce$value_order), ]
trust_amce$value <- factor(trust_amce$value, levels = rev(unique(as.character(trust_amce$value))))

## Plot the AMCEs
fig_trust <- ggplot(trust_amce, aes(x = estimate, y = value)) +
  facet_grid(treatment ~ .,
             scales = "free_y",
             space = "free_y") +
  geom_errorbarh(aes(xmin = estimate - (2 * std_error),
                     xmax = estimate + (2 * std_error)),
                 height = 0) +
  geom_point() +
  geom_vline(aes(xintercept = 0),
             linetype = "dotted") +
  scale_x_continuous(limits = c(-.3, .3),
                     breaks = round(seq(-.3, .3, .1), 1),
                     expand = c(0, 0),
                     labels = function(x) x * 100) +
  scale_y_discrete(labels = function(x) parse(text = as.character(x))) +
  labs(x = "Marginal effect, choosing source (%)",
       y = NULL) +
  theme_m() +
  theme(strip.text.y = element_blank())

ggsave("fig_trust.pdf", plot = fig_trust, width = 7, height = 6)



## Example 2: Selective exposure conjoint ----------------------------

## Get the data.
selective <- read.csv("data-selective.csv", stringsAsFactors = FALSE)

## Fix spelling errors.
selective$hl_opinion <- gsub("Religous", "Religious", selective$hl_opinion)
selective$hl_opinion <- gsub("servies", "services", selective$hl_opinion)

## Labels to use for treatment headers in the figure.
selective_labels <- data.frame(
  treatment = c(
    "hl_opinion",
    "hl_direction",
    "hl_party_nor",
    "hl_matched_opinion",
    "hl_matched_party_value_cat"
  ),
  treatment_label = c(
    "Message topic",
    "Message direction",
    "Party mentioned",
    "Agree with message",
    "Likes party mentioned"
  ),
  stringsAsFactors = FALSE
)


## Make the vars factors.
for (var in selective_labels$treatment)
  selective[[var]] <- factor(selective[[var]])

## For the purpose of this paper, remove the "Neither" condition for matched message and party.
selective$hl_matched_opinion <- ifelse(selective$hl_matched_opinion == "Neither", NA, as.character(selective$hl_matched_opinion))
selective$hl_matched_party_value_cat <- ifelse(selective$hl_matched_party_value_cat == "Neither", NA, as.character(selective$hl_matched_party_value_cat))

## Fix treatment levels and order as per wanted in the figure.
selective$hl_opinion <- factor(selective$hl_opinion, levels = levels(selective$hl_opinion)[c(3, 2, 1, 5, 6, 4)])
selective$hl_party_nor <- factor(selective$hl_party_nor, levels = levels(selective$hl_party_nor)[c(6, 8, 1, 7, 5, 9, 4, 3, 2)])
selective$hl_matched_opinion <- factor(selective$hl_matched_opinion, levels = unique(selective$hl_matched_opinion)[c(2, 1)])

## Estimate the AMCEs.
selective_amce <- amce(
  data = selective,
  post = "hl_choice",
  treatments = c(
    "hl_opinion",
    "hl_direction",
    "hl_party_nor",
    "hl_matched_opinion",
    "hl_matched_party_value_cat"
  ),
  cluster = "resp_id"
)

## Adding a row with treatment labels so the plot have a header on the y axis.
selective_amce <- selective_amce[rev(order(selective_amce$value)), ]
selective_amce$value <- paste0("plain('", as.character(selective_amce$value), "')")
selective_amce$value <- factor(selective_amce$value, levels = unique(selective_amce$value))
add_data <- selective_labels
names(add_data)[2] <- "value"
add_data$value <- gsub("\n", " ", as.character(add_data$value))
add_data$value <- paste0("bold('", add_data$value, "')")
for (var in c("value_order", "estimate", "std_error", "statistic", "p_value")) {
  add_data[[var]] <- NA
}
add_data$value_order <- -Inf
selective_amce <- rbind(selective_amce, add_data)
selective_amce <- selective_amce[order(selective_amce$treatment, selective_amce$value_order), ]
selective_amce$value <- factor(selective_amce$value, levels = rev(unique(as.character(selective_amce$value))))

## ## Plot the main AMCEs
## fig_selective <- ggplot(selective_amce, aes(x = estimate, y = value)) +
##   facet_grid(treatment ~ .,
##              scales = "free_y",
##              space = "free_y") +
##   geom_errorbarh(aes(xmin = estimate - (2 * std_error),
##                      xmax = estimate + (2 * std_error)),
##                  height = 0) +
##   geom_point() +
##   geom_vline(aes(xintercept = 0),
##              linetype = "dotted") +
##   scale_x_continuous(limits = c(-.25, .25),
##                      breaks = round(seq(-.25, .25, .05), 2),
##                      expand = c(0, 0),
##                      labels = function(x) x * 100) +
##   scale_y_discrete(labels = function(x) parse(text = as.character(x))) +
##   labs(x = "Marginal effect, choosing headline (%)",
##        y = "Headline attributes") +
##   theme_m() +
##   theme(plot.margin = unit(c(2, 2, 2, 2), "mm"))
## fig_selective

## Plot the AMCEs of the matched attributes.

## Plot the main AMCEs
fig_selective <- ggplot(selective_amce[1:20, ], aes(x = estimate, y = value)) +
  facet_grid(treatment ~ .,
             scales = "free_y",
             space = "free_y") +
  geom_errorbarh(aes(xmin = estimate - (2 * std_error),
                     xmax = estimate + (2 * std_error)),
                 height = 0) +
  geom_point() +
  geom_vline(aes(xintercept = 0),
             linetype = "dotted") +
  scale_x_continuous(limits = c(-.25, .25),
                     breaks = round(seq(-.25, .25, .05), 2),
                     expand = c(0, 0),
                     labels = function(x) x * 100) +
  scale_y_discrete(labels = function(x) parse(text = as.character(x))) +
  labs(x = "Marginal effect, choosing headline (%)",
       y = "Headline attributes",
       title = "a) Randomized headline attributes\n") +
  theme_m() +
  theme(plot.margin = unit(c(2, 2, 2, 22), "mm"),
        plot.title = element_text(size = 12, face = "italic", hjust = 0))
fig_selective

## Plot the main AMCEs
fig_selective_matched <- ggplot(selective_amce[21:26, ], aes(x = estimate, y = value)) +
  facet_grid(treatment ~ .,
             scales = "free_y",
             space = "free_y") +
  geom_errorbarh(aes(xmin = estimate - (2 * std_error),
                     xmax = estimate + (2 * std_error)),
                 height = 0) +
  geom_point() +
  geom_vline(aes(xintercept = 0),
             linetype = "dotted") +
  scale_x_continuous(limits = c(-.25, .25),
                     breaks = round(seq(-.25, .25, .05), 2),
                     expand = c(0, 0),
                     labels = function(x) x * 100) +
  scale_y_discrete(labels = function(x) parse(text = as.character(x))) +
  labs(x = "Marginal effect, choosing headline (%)",
       y = "Headline attributes\nmatched with\nrespondent",
       title = "b) When matched with respondent preferences\n") +
  theme_m() +
  theme(plot.margin = unit(c(10, 2, 2, 26.125), "mm"),
        plot.title = element_text(size = 12, face = "italic", hjust = 0))
fig_selective_matched

## Print to figure to file.
pdf("fig_selective.pdf", width = 7, height = 6.5)

grid.arrange(fig_selective, fig_selective_matched,
             heights = c(6.5/10, 3.5/10),
             ncol = 1)

dev.off()


## Numbers in text ---------------------------------------------------

## ## TRUST: Number if respondents
## length(unique(trust$resp_id))

## ## TRUST: Number of choice tasks
## nrow(trust)

## ## SELECTIVE: Number of unique headlines
## length(unique(selective$hl_headline))

## ## SELECTIVE: Number of resp
## length(unique(selective$resp_id))

## ## SELECTIVE: Number of choice task
## nrow(selective)
